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PART  It  The  6-Layer  AWSPE  Model 


Major  Terry  C.  Tarbell 
Captain  Fred  P.  Lewis 


1.  INTRODUCTION 

The  Air  Weather  Service  Primitive-Equation  (AWSPE)  model,  also 
referred  to  as  the  6-Layer  (6L)  AWSPE,  is  a  dry  version  (that  is,  no  moisture 
variable  is  forecast)  of  the  Shuman-Hovermale  model  (Shuman  and  Hovermale, 
1968).  Before  discussing  the  history  of  the  6L  AWSPE,  it  is  appropriate  to 
present  a  brief  history  of  the  Shuman-Hovermale  model. 

1.1  The  Shuman-Hovermale  Model  at  NMC. 

On  6  June  1966,  the  Shuman-Hovermale  six-layer  primitive-equation  (6L 
PE)  model  became  operational  at  NMC  (this  entire  subsection  was  adapted  from 
Shuman  and  Hovermale,  1968) .  The  6L  PE  was  the  first  primitive-equation  (PE) 
model  used  operationally  anywhere. 

The  6L  PE  model  replaced  the  Cressman  (1963)  three-level  filtered 
equation  model.  The  Cressman  model  and  its  predecessors  at  NMC  were  all  based 
on  Rossby's  (1939)  idea  that  the  principle  mechanism  in  the  short-term 
prediction  of  large-scale  atmospheric  motions  is  the  conservation  of  the 
vertical  component  of  absolute  vorticity.  The  application  of  Rossby's  idea 
leads  to  the  use  of  differentiated  forms  of  the  equations  of  motion  rather 
than  the  "primitive”  equations  themselves.  Charney  (1962)  neglected  the 
time-rate-of-change-of-divergence  term  in  the  divergence  equation.  The 
resulting  equations  are  called  "balance  equations”  and  these  equations 
"filter"  or  eliminate  the  fastest  moving  waves  (inertia-gravity  waves) 
permitted  by  the  primitive  equations.  Thus,  the  only  motions  permitted  by  the 
filtered  equations  are  the  slowly  moving  meteorological  waves  of  interest. 

By  1959,  NMC  scientists  realised  that  further  advances  in  numerical 
weather  prediction  (NWP)  would  be  severely  limited  unless  the  primitive 
equations  themselves  were  used.  This  led  to  the  development  of  the  6L  PE 
model  between  1959  and  1966.  Many  obstacles,  such  as  computational 
instability,  had  to  be  overcome  in  the  development  of  the  model.  After  the  6L 
PE  was  developed,  it  did  not  meet  operational  timelines  because  it  required 
many  more  calculations  than  the  Cressman  model.  However,  two  events  occurred 
that  made  operational  use  of  the  model  feasible: 

a.  The  development  of  a  new  generation  of  computers  and  the  installation 
of  one  of  them  (a  CDC-6600)  at  NMC. 

b.  The  Automated  Weather  Network  (AWN)  that  collected  data  significantly 
faster  than  before,  thereby  permitting  an  earlier  6L  PE  start  time.  The 
operational  implementation  of  the  Shuman-Hovermale  model  will  always  be 
remembered  as  a  major  achievement  in  the  field  of  numerical  weather  prediction 


1.2  Performance  of  the  Shuman-Hovermale  Model 


Shuman  and  Hovermale  (1968)  presented  the  S-l  scores  (Teweles  and 
Wobus,  1954)  for  the  months  between  June  1966  and  August  1967  for  the 
barotropic  model  36-h  500-mb  forecasts  and  for  the  61.  PE  36-h  500-mb 
forecasts.  The  6L  PE  model  set  the  record  low  S-l  score  for  any  month  four 
times  and  set  records  for  seven  months  of  the  year.  The  6L  PE  S-l  scores  beat 
the  barotropic  model  S-l  scores  by  an  average  of  6  points  for  the  15-month 
verification  period.  Since  an  improvement  in  SI  score  of  4  points  is 
considered  good,  an  improvement  of  6  points  is  indeed  impressive. 

Many  other  statistics  also  led  to  the  conclusion  that  the  6L  PE  was 
far  superior  to  the  Cressman  model.  Thus,  the  use  of  primitive-equation 
models  to  forecast  the  atmosphere  was  firmly  established. 

1.3  A  6-Layer  Primitive-Equation  Model  at  AFGWC. 

The  6-Layer  AWSPE  (6L  AWSPE) ,  a  version  of  the  Shuman-Hovermale  model, 
was  implemented  at  APGWC  on  4  May  1975.  The  6L  AWSPE  was  the  first 
primitive-equation  model  to  be  used  at  APGWC.  It  currently  runs  four  times 
daily  in  the  Northern  Hemisphere  producing  72-h  forecasts  for  the  00  and  12 
GMT  data  times  and  36-h  forecasts  for  the  06  and  18  GMT  data  times.  Tarbell 
and  Hoke  (1979)  presented  a  description  of  how  the  6L  AWSPE  fits  in  the 
operational  production  system  at  APGWC. 

The  previous  operational  forecast  model  at  APGWC  was  SIXLVL,  the 
Six-Level  Quasi-geostrophic  Forecast  Model,  described  by  Palucci  (1970) . 

Unlike  the  6L  AWSPE,  SIXLVL  is  a  filtered  model.  Flattery  (1975)  reported 
that  the  formulation  of  SIXLVL  is  based  on  a  number  of  assumptions  and 
approximations  including: 

a.  The  horizontal  wind  is  assumed  to  be  nondivergent. 

b.  The  static  stability  is  a  function  of  pressure  only. 

c.  An  artificial  term  is  included  to  stabilize  the  long  waves. 

d.  Surface  friction  is  included  only  in  the  form  of  a  modified  vertical 
velocity. 

e.  The  Coriolis  parameter,  f,  is  assumed  to  be  constant  for  certain 
equation  terms. 

f.  Some  smaller-magnitude  terms  have  been  omitted  from  the  model 
equations. 

The  6L  AWSPE  model  does  not  use  any  of  these  assumptions.  Thus,  after 
implementing  the  6L  AWSPE  model,  APGWC  realized  improved  forecasts  similar  to 
the  improved  forecasts  that  occurred  at  NMC  after  implementation  of  the  6L 
PE.  Major  differences  between  6L  AWSPE  and  the  SIXLVL  combination  can  be 
expected  in  vertical  velocity  and  vorticity.  In  the  SIXLVL  model,  vertical 
velocity  fields  are  dependent  only  on  the  stream  function.  In  the  6L  AWSPE 
model,  the  vertical  velocity  field  is  derived  directly  from  the  total  rate  of 
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change  with  respect  to  time  in  the  pressure  field  and  therefore  can  be  quite 
different  in  both  detail  and  Magnitude. 


In  the  SIXLVL  model,  the  winds  are  diagnosed  and  remain  nearly 
geostrophic.  In  the  6L  AWSPE  model,  the  wind  components  themselves  are 
forecast  and  significant  ageostrophic  components  result.  Thus,  the  6L  AWSPE 
vorticity  centers  can  be  displaced  from  the  positions  that  would  be  expected 
from  the  height-contour  pattern. 

Since  SIXLVL1 8  lowest  forecast  level  is  850  mb,  another  forecast  model  is 
used  to  provide  surface  forecasts.  Originally  developed  by  Reed  (1963) , 

1MPROG  produces  the  1000-mb  forecast  by  advecting  the  1000-500  mb  thickness 
with  a  portion  of  the  500-mb  wind  field.  In  the  6L  AWSPE  model,  surface 
forecasts  are  obtained  through  the  dynamic  evolution  of  the  three-dimensional 
model  atmosphere  using  a  terrain-following  coordinate  system. 

The  remainder  of  Part  I  of  this  Tech  Note  gives  a  technical  description  of 
the  6L  AWSPE  model.  Section  2  contains  a  description  of  the  vertical  and 
horizontal  model  domains  and  the  set  of  primitive  equations  used.  In  Section 
3,  we  present  the  model  finite-difference  equations  and  other  physical 
parameter izations.  In  Section  4,  the  sequence  of  calculations  used  to  compute 
the  forecast  are  detailed.  Section  5  contains  the  procedure  used  to 
initialize  the  6L  AWSPE.  In  Section  6,  we  present  the  procedure  used  to 
calculate  the  forecast  fields  on  mandatory  pressure  surfaces.  Section  7 
contains  a  comparison  of  the  6L  AWSPE  with  the  original  Shuman-Hove rmale 
model.  Section  8  describes  the  subroutines  in  the  6L  AWSPE  and  contains  a 
basic  subroutine  flow  diagram  for  the  model.  This  section  is  mainly  for  use 
by  AFGWC/TSIN  programmers  and  may  be  skipped  by  other  readers.  Section  9  is  a 
summary  of  Part  I. 

2.  ASWPE  MODEL  CHARACTERISTICS 


The  Primitive 


lations. 


The  hydrostatic  meteorological  equations  for  a  model  with  a  vertical 
sigma  coordinate  and  horizontal  Cartesian  coordinates  on  a  conformal  polar 
stereographic  projection  of  the  earth  are  given  by  Shuman  and  Hovermale 
(1968) .  They  are: 
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where  o  is  defined  ui 


a  -  P  ~  Py 
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(2.7) 


These  six  equations  ere  written  in  the  six  dependent  variables  u,  v,  o, 

9  ,  it  ,  and  p.  It  is  worth  noting  that  x  and  y  are  the  two  horisontal 
Cartesian  coordinates  on  the  projection  while  u  and  v  are  the  x  and  y 
ccsqponents  of  velocity  true  on  the  earth.  The  sap  factor,  a,  is  the  ratio  of 
a  distance  on  the  projection  to  the  corresponding  distance  on  the  earth. 
Therefore,  we  can  write 


°proj 

vproj 


(2.8) 
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wiwrt  OpCOj  and  Vproj  a ra  tha  horizontal  velocity  components  cm  the 
conformal  projection. 

The  definition  of  sigma  given  by  (2.7)  ie  a  generalisation  of  Phillips 
(1957)  vertical  coordinate  e  ■  p/p*,  where  p*  ia  surface  pressure. 

2.2  The  vertical  Domain. 

Pig.  2.1  is  a  depiction  of  the  vertical  structure  of  the  6L  AMSPE 
model.  The  model  contains  four  separate  sigma  dooming  that  contain  the  seven 
vertical  layers.  The  lowest  sigma  domain  is  a  planetary  boundary  layer  with 
the  vertical  coordinate  defined  by 
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PsPs 


The  second  sigma  domain  consists  of  three  layers  and  constitutes  the 
remainder  of  the  troposphere.  The  second  sigma  doawin  is  defined  by 


Oj  -  P  ~ 

P5  -  P? 


(2. 11^ 


where  p£*  is  the  pressure  at  the  tropopause.  The  tropopauae  is  considered  to 
be  a  material  surface;  that  ia,  *  °* 

The  third  sigma  domain  consists  of  two  layers  and  constitutes  the  model 
stratosphere.  The  stratospheric  domain  is  defined  by 


p\*-Po  (2,12> 

where  pQ  is  the  pressure  at  the  top  of  the  model  (*100mb) . 

The  level  at  the  top  of  the  stratosphere  (k-0)  is  also  treated  as  a 
material  surface. 

The  topmost  sigma  domain  is  a  layer  of  constant  potential  teaperature.  it 
is  defined  by 


o  -  _E 
e  p9 


(2.13) 


This  seventh  lsyer  is  included  for  its  coaputational  rather  than  its 
meteorological  properties. 

2.3  The  Horisontal  Grid. 

The  horisontal  grid  is  a  regular  square  aesh  on  a  polar  stereographic 
projection  true  at  SON.  For  this  projection. 
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(2.14) 


where  0  is  north  latitude.  Centered  at  the  North  Pole,  the  horisontal  grid 
consists  of  S3  X  57  (3021)  points  (see  Fig.  2.2).  The  grid  increaent  is  one 
bedient  (381  ka)  at  60N.  The  APGNC  Northern  Beaisphere  octagon  grid  is  a 
subset  of  points  within  the  ANSPB  horisontal  grid  (see  Hoke  at  al. .  1981). 
Model  forecasts  produced  outside  the  octagon  are  not  used.  The  boundary 
conditions  at  the  edge  of  the  rectangular  grid  and  the  initialisation 
technique  eaployed  between  the  octagon  and  the  rectangular  grid  caabine  to 
produce  less  reliable  forecasts  outside  the  octagon.  The  lateral  boundary 
conditions  and  the  initialisation  procedure  are  dlacussed  in  later  sections. 
The  purpose  of  the  boundary  region  is  to  provide  a  buffer  sone  between  the 
walls  and  the  active  Meteorological  areas. 

2.4  Prognostic  and  Diagnostic  Variables. 

The  forecast  or  history  variables  are  u,  v,  9,  and  p  .  The  diagnosed 
variables  are  if,  i,  6,  and  p.  In  Fig.  2.3,  we  show  that  u?  v,  and  a  are 
defined  at  the  half  levels  and  p,  ir,  6,  and  s  are  defined  at  the  full  levels. 
Actually,  o  is  defined  at  the  center  point  of  each  square  foraed  by  four  grid 
points  at  which  s,  ir ,  and  p  are  defined  (see  Fig.  2.3) .  The  Methods  used  to 
calculate  these  variables  are  presented  in  sections  3  and  4. 

3.  FINITE-DIFFERENCE  FORMULATION  OF  TEE  ANSPB  MODEL. 

In  this  section,  we  will  present  the  aodel  finite-difference 
equations,  the  PBL  friction  paraaeterisation,  the  aethod  of  calculating  a, 
the  lateral  boundary  conditions,  the  horisontal  saoothing  and  the  tine 
snoot hing,  the  pressure-gradient-averaging  Methods,  and  the  techniques  used  to 
coapute  oaega  (ui),  which  is  the  vertical  notion  in  pressure  coordinates. 


ANSPB  Finite-Difference 


tions. 
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Figure  2.2. 


The  AFGWC  northern  hemisphere  whole-mesh  grid  for  the 
6L  and  7L  AWSPE  models  Is  given  by  the  Inner  solid 
lines.  The  domains  of  the  Northern  Hemispheric 
Reference  Grids  (solid  outer  border)  and  the  Octagon 
rids  (dashed  lines)  are  also  Indicated.  The  Indices 
I  ,1  )  designate  the  coordinates  for  both  the 
Whole-mesh  Reference  Grid  and  the  AWSPE  grid.  The 
whole  mesh  grid  spacing  Is  displayed  In  the  upper- 
left  comer. 
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where  i  is  tbs  ssrisl  number  of  «  point  in  ths  x  dirsetlon.  The  use  of  more 
than  one  superscript  or  subscript  neans 


(3.3) 
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(3.4) 


Values  of  f  are  defined  at  grid  points  on  the  horisontal  grid.  Therefore, 
values  of  fx  and  f*  represent  values  between  grid  points.  The  following 
combinations  of  operations  are  axaaples  of  or— on  finite-difference  operations t 
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Using  the  ShuMn-Boveraale  finite-difference  notation,  the 
finite-difference  representations  of  (2.1)  through  (2.5)  are  written t 
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Tht  (  )t  ovir  th«  pressure-gradient-force  Urn  In  (3.7)  and  (3.8) 
represents  the  pr«MHir«-amaging  technique  of  Brown  and  Caapana  (1978) . 

In  this  technique,  tha  pressure-gradient  forca  ia  coaputed  at  thraa  tiaa  steps 
(t+1,  t,  and  t-1,  wbara  t  la  tha  tiaa  atap)  and  than  averaged  aa  follows: 


)*  -  uK  )*+l  ♦  (  I*"1  ]  ♦  (?  -  fa)  (  )*.  (3.12) 


This  tachniqua  is  nuaerleally  atatola  for  a  ■  0.25.  Xn  practice,  a  »  0.245 
for  tha  first  six  for scant  hours  and  is  than  sat  at  a*  0.25  for  tha  reaainder 
of  tha  foracast.  Tha  pressure-gradient-averaging  tachniqua  allows  tha  tiaa 
atap  used  by  tha  nodal  to  ba  doublad. 


3.2  PBL  Friction. 


Frictional  dissipation  is  only  considered  in  tha  aodel  boundary 
layer.  Thus,  vortical  aonontua  exchange  due  to  turbulent  (subgrid-scale) 
addles  is  paraaeteriied  in  tha  nodal  boundary  layer.  The  friction  terns  rx 
and  Fy  in  (3.7)  and  (3.8)  have  the  following  forast 
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Hhttt  p  is  tbs  standard  sea- level  danaity  and  dp  ■  SO  ab.  lbs  drag 
coefficient,  Cq,  is  aodslsd  aftsr  Cressaan  (I960)  and  has  a  Xarga  horizontal 
variation  that  is  positivoly  corrslatad  with  tarrain  haight. 

3.3  Tha  Calculation  of  "Sigaa  Dot*. 

Tha  aathod  usad  Cor  coaputing  «  follows  alaoat  diractly  Croa  tha 
finite. difference  fora  of  tha  p0  tandancy  aquation  (3.11)  and  tha  fact 
that  a  -  0  at  lavals  a  -  6*#  2**$  0 ,  and  -1. 


First,  p »  50  ab  in  tha  boundary  layar  and  is  constant  in  tiaa  and 
spaca.  Thus,  if  wa  raaova  tha  PT***  operator,  (3.11)  nay  be  written  as 
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where 
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in  tha  troposphere,  a  different  approaeh  auat  be  usad,  since  p0  varies  in 
tiaa  and  spaca.  However,  p0  auat  be  linear  with  respect  to  o  ,  therefore 
Poo  ■  °*in9  this  fact,  an  aquation  for  So0  °*n  be  found  and  than 

(dT) 4  and  (dT)3  can  be  ccaputad.  Differentiating  (3.11)  with  respect  to 
yields 
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then  (3.19)  can  be  written  es 
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Applying  (3.20)  at  levels  k  ■  3  and  4  yields 
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Solving  (3.21)  for  (0.^3  and  (oT>4  given  that 
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Siailarly,  in  the  stratosphere 


3.4  Lateral  Boundary  Condition. 

The  6L  ANSPE  has  rigid  impermeable  walls  between  the  two  outermost 
grid  rows  (see  Fig.  2.2).  No  heat  or  momentum  is  allowed  to  pass  through  the 
boundary,  although  flow  parallel  to  the  boundary  is  permitted.  In  finite- 
difference  fora  these  conditions  may  be  written  as: 


for  x-constant 
and  for  y  *  constant 

U,.  *  0  i  e„  *  0  ,  \>V  •  0  ,?•  0  .  (3.28) 

y  '  y  '  °y 


These  conditions  lead  to  the  reflection  of  pure  gravity  waves  at  the 
boundaries. 

3.5  Theraodynaaic  Forcing. 

The  eddy  flux  of  heat  froa  ocean  surfaces  is  parameterized  after  the 
work  of  Lorenz  (1962).  The  boundary  layer  is  only  heated  if  the  sea  surface 
temperature  is  greater  than  the  boundary-layer  potential  temperature 
(8).  The  6L  ANSPE  model  uses  mean  monthly  sea  surface  temperature  fields  to 
derive  the  heating  rate.  The  boundary-layer  heating  rate  is  given  by 

H  *  K  (8  -  TlXy 

w  (3.29) 

where  K  «  10“4  s-1  and  Tw  is  the  sea  surface  temperature. 

3.6  Horizontal  Smoothing. 

Explicit  horizontal  space  smoothers  are  not  used  in  the  central 
portion  of  the  model  domain  (north  of  20N) .  Below  20N,  two  types  of  smoothers 
(tendency  truncation  and  lateral  diffusion)  are  used. 

The  tendency  truncation  is  applied  to  the  u,  v,  8,  and  pa  equations.  This 
form  of  smoothing  is  used  to  reduce  noise  in  the  Tropical  regions  of  the 
model.  The  model  uses  half  of  the  computed  tendency  for  the  above  variables 
below  ION  and  the  full  tendency  at  and  above  20N.  The  two  zones  are  blended 
linearly  between  ION  and  20N. 

The  lateral  diffusion  "smoother"  is  introduced  into  the  u  and  v  momentum 
equations  as  follows: 

~  ■  (u  tendency  terms)  +  A‘ (KAu)  (3.30) 
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3.7  Time  Smoothing. 

Because  the  model  uses  centered  (leapfrog)  tine  differencing,  a  time 
smoother  is  applied  to  control  time  splitting.  The  time  smoother  has  the  form: 


24>n  +  tiF1*1  ) 


(3.33) 


where  4*  represents  the  variables  u,  v,  0,  and  p0;  n  represents  a  given  time 
level;  and  a  «*  0.075  for  even  time  steps  and  -0.075  for  odd  time  steps. 

3.8  Computation  of  Omega. 


Omega  (ui)  is  computed  only  for  output  for  use  by  forecasters  and 
other  models.  The  model's  vertical  motion  parameter  is  "sigma  dot"  (the  total 
time  rate  of  change  of  sigma)  since  the  AMSPB  model  was  formulated  using  sigma 
coordinates. 

Omega  is  computed  directly  using 


or  in  finite  difference  form 
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(3.35) 


The  value  of  w  is  then  averaged  over  2  h  (6  time  steps  for  which  at  ■  1200  s) 
before  output.  For  instance,  the  omega  value  "valid"  at  the  3-h  forecast 
point  would  be  an  average  of  the  values  between  forecast  hour  1  and  forecast 
hour  3. 


Computation  of  the  Geopotential  Height 


The  geopotential  height  £or  the  6L  AWSPE  is  computed  directly  from 
the  finite-difference  formulation  of  the  hydrostatic  approximation.  Eg.  (3.9) 
is  integrated  in  the  vertical,  given  the  terrain  height,  potential 
temperature,  and  pressure  of  the  sigma  levels.  Values  of  $  computed  at  the 
whole  levels  are  averaged  to  the  half  levels  where  they  are  used  in  computing 
the  u  and  v  tendencies. 

3.10.  Dry  Convective  Adjustment. 

The  6L  AWSPE  model  applies  a  dry  convective  adjustment  procedure  to 
remove  superadiabatic  lapse  rates  and  thus  maintain  model  stability.  Haltiner 
and  williams  (1980)  show  that  the  phase  speed  of  internal  gravity  waves  is 
given  by 
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where  m  and  k  are  the  wave  numbers  in  x  and  z,  respectively,  g  is  gravity, 
and  6  is  the  mean  value  of  the  potential  temperature  using  a  standard 
perturbation  analysis  with 


e  *  8 (z)  +  e’(x,z,;t) 


(3.37) 
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When  the  vertical  lapse  rate  is  superadiabatic,  3z  becomes  negative  and  the 
phase  speed  is  complex.  Thus,  the  internal  gravity  wave  is  unstable  and  will 
grow  without  bound;  i.e.  the  model  will  become  computationally  unstable. 
Therefore,  dry  convective  adjustment  is  applied  in  the  model  to  adjust  all 
superadiabatic  lapse  rates  to  the  adiabatic  or  neutral  position  with  39  ■  0. 

3z 

4.  THE  6-LAYER  AWSPE  MARCHING  PROCESS. 

The  finite-difference  equations  are  solved  as  follows: 

a.  Given  initial  values  for  u,  v,  e,  and  po  ,the  future  values  of 
these  variables  are  computed  using  (3.7)  through  (3.11).  In  order  to  apply 
the  pressure-gradient-averaging  technique,  the  future  values  of  pa  and  e  must 
be  computed  before  the  pressure-gradient  force  can  be  calculated.  First  o  is 
computed  using  the  method  given  in  Section  3.3.  p  is  now  forecast  from 
(3.11)  and  0  is  computed  using  (3.10).  Note  that  0  is  not  forecast  in  the 
isentropic  top  layer.  Next,  the  pressure  valuer  ■»':  levels  0,  1,  2,  3,  4,  5, 


and  6  for  the  future  tine  step  N+l  (where  N  represents  an  arbitrary  tine  level 
or  step)  are  computed  by  summing  p  downward  from  the  top  level  (-1)  where  p 
-  0. 

b.  Exner's  function,  it  ,  at  time  step  N+l  is  computed  from 
(2.6).  Then  (3.9)  is  used  to  compute  the  height  fields,  z,  at  levels  5 
through  -1  for  time  N+l  given  the  terrain  field  at  level  6.  The 
pressure-gradient  force  can  be  computed  from  the  values  of  z  ,  0,  and  it  at  time 
steps  N-l,  N,  N+l  using  (3.7),  (3.8),  and  (3.12). 

c.  Eqs.  (3.7)  and  (3.8)  are  then  solved  for  u  and  v  at  time  step 

N+l. 

This  process  is  repeated  in  time  until  the  desired  forecast  time  is  reached. 
For  the  first  time  step,  the  model  uses  forward  differencing.  After  the  first 
time  step,  a  centered  difference  scheme  is  applied.  Because  the  model  uses 
pressure-gradient  averaging  a  time  step  of  1200s  is  permitted. 

5.  INITIALIZATION 

The  initialization  process  is  detailed  in  Shuman  and  Hovermale 
(1968).  Currently,  heights  on  pressure  coordinates,  taken  from  the  APGMC 
Hough  Spectral  Analysis  Model  (HUFANL) ,  are  interpolated  to  the  6L  AMSPB 
sigma-coordinate  domain  using  the  temperature  fields  from  the  Hough  Spectral 
Analysis  Model  (Tarbell  and  Hoke,  1979).  These  height  values  are  then  used  to 
compute  potential  temperature,  6  .  The  procedure  uses  a  linear  balance 
equation  to  compute  the  winds  on  the  sigma  domain  given  the  height  and 
potential  temperature  fields. 

This  balanced  initialization  procedure  typically  reduces  the  mass-weighted 
kinetic  energy  by  about  10  to  15  percent  when  compared  to  the  mass-weighted 
kinetic  energy  computed  from  the  actual  analysis. 

6.  OUTPUT  FIELDS 

The  model  forecast  values  on  sigma  surfaces  are  stored  in  a  file  on 
the  computer  every  3h  of  forecast  time.  A  post  processor  interpolates  the  3h 
sigma  data  on  sigma  surfaces  to  the  AFGMC  pressure-coordinate  data  base.  This 
interpolation  procedure  computes  the  value  of  u,  v,  and  T  for  the  mandatory 
pressure  levels  by  linear  interpolation  of  the  sigma- level  values  in  Exner's 
function  (it  ).  The  height  field  (z)  at  the  sigma  levels  is  computed  from 
the  hydrostatic  equation.  Then,  the  height  and  temperature  fields  are 
interpolated  in  it  to  mandatory  pressure  levels.  The  exact  details  of  this 
procedure  can  be  found  in  the  6L  AWSPE  Program  Maintenance  Manual  (Andrew 
Johnson,  1975) . 

The  pressure-coordinate  data  base  is  then  stored  in  the  AFGHC 
literal-label  data  base  (under  literal  RT6MCF) .  This  data  base  is  used  by 
display  programs  and  other  models. 


This  literal-label  data  bass  contains  wind,  taaparaturs,  and  D-value 
fields  in  3-  to  6-h  increaents  frea  ths  surface  to  100  ab.  Fields  such  as 
surface  pressure,  vorticity,  vertical  notion  (oaega) ,  and  tropopause  height 
and  teaperature  are  carried  at  certain  levels  and  tines.  Values  are  routinely 
forecast  in  the  Northern  Reaiisphere  to  72h  for  the  00  and  12  GMT  data  tines 
and  to  3<h  for  the  06  and  18  GMT  data  tines.  Table  7.1  (froa  the  AFGMC  Data 
Base  Handbook)  lists  the  fields  stored  in  the  6L  AM8PB  data  base. 


18 


TABLE  6.1.  6L  AMS PE  Data  Basa 


FIELD 


Label  Laval  (LL) tab* 10)  Forecast  Tiae  <TT) (h) 


D-value 

TTAXLL 

00 

03, 

06, 

•  •  •  9 

24, 

30, 

36, 

48, 

60, 

72 

TTAXLL 

IS. 

70, 

50, 

40,  30, 

06, 

12, 

•  •  •  9 

48, 

60, 

72 

25, 

20, 

15, 

10 

TTAXLL 

50 

96 

Wind 

TTAMLL 

00 

00, 

03, 

06, 

•  •  •  9 

60 

TTAWLL 

•5, 

70, 

50. 

30, 

00, 

03, 

•  ••# 

60, 

72 

20, 

10 

TTAMLL 

60, 

25, 

15 

06, 

12, 

•  ••# 

48, 

60, 

72 

TTAMLL 

OR 

09 

Sfc  Pres/ 

TTAPLL 

BF 

03, 

06, 

•  •  •  9 

24, 

30, 

36 

lOOOab  reap 

Vertical  Velocity 

TTADLL 

70, 

30, 

10 

00, 

03, 

•  •  •  9 

60 

Streaa  Function 

TTASLL 

05, 

70, 

50, 

00, 

03, 

•  •  •  9 

48 

30, 

20, 

10 

Tropopause 

TTALLL 

Tropopause 

06, 

12, 

•  •  •  9 

48, 

60, 

Pres/Hgt/Tenp 

72 

Vorticity 

TTAVLL 

50 

06, 

12, 

•  •  •  9 

48 

7.  DIFFERENCES  BETWEEN  THE  6-LAYER  AWSPE  MODEL  AND  THE  IMC  7-LAYER 
SHUMAN- HOVERMALE  PE  MODEL 

The  6L  AWSPE  is  an  early  version  of  the  7-Layer  Shuman-Hovermale  pe 
nod el  (7  LPE)  used  operationally  at  NMC  until  August  1980.  The  basic 
structure  of  the  two  models  is  very  similar ,  but  substantial  differences  do 
exist.  In  fact,  the  7L  PE  forecasts  verify  substantially  better  than  the  6L 
AWSPE  forecasts. 

Fig.  7.1  shows  root-mean-square  vector  errors  at  500  ab  for  the  two 
models,  the  1-bedient  5L  AWSPE  and  the  1/2-bedient  7L  PE  for  the  period 
January  1979  to  December  1980.  One  bedient  is  exactly  equal  to  what  we  at 
AFGWC  refer  to  as  whole  mesh,  see  Hoke,  et.al.  (1980)  for  further 
explanation.  This  figure  clearly  demonstrates  the  quality  of  the  enhancements 
made  by  NMC  to  their  version  of  the  1966  Shuman-Hove rmale  model. 

We  might  ask  why  does  the  1/2-bedient  7L  PE  verify  better  than  the  6L 
AWSPE.  Haltiner  and  Williams  (1980)  pointed  out  that  reduced  phase-speed 
error  in  the  weather-producing  waves  leads  to  improved  short-range  forecasts 
for  models  like  the  7L  PE.  They  stated  phase  speed  errors  can  be  reduced  by 
minimising  truncation  errors.  Thus,  one  of  the  most  important  reasons  why  the 
6L  AWSPE  doesn't  produce  forecasts  as  good  as  the  NMC  7L  PE  model  is  probably 
the  reduced  phase-speed  errors  due  to  the  1/2-bedient  grid  used  by  the  NMC  7L 
PE.  Less  important  reasons  exist.  For  example,  the  7L  PE  is  a  moist  model 
that  parameterises  latent  heating  effects.  The  7L  PE  also  uses  a  more 
sophisticated  initialisation  technique. 

8.  A  BRIEF  DESCRIPTION  OF  THE  SUBROUTINES  IN  THE  6-LAYER  AWSPE  PROGRAM 

In  this  section,  we  will  present  a  brief  description  of  the  different 
subroutines  in  the  6L  AWSPE  model  taken  from  Johnson  (1975).  This  section  is 
intended  for  those  requiring  a  more  detailed  knowledge  of  the  6L  AWSPE. 

Pig.  8.1  gives  the  6-Layer  AWSPE  model  subroutine  flow  diagram.  The 
numbers  in  the  lower-left  hand  corner  of  each  subroutine  give  the  order  in 
which  the  programs  are  used.  The  basic  function  of  each  module  shorn  in  the 
diagram  is  defined  below. 

UTIME3  -  The  main  routine  for  the  6L  AWSPE  forecast  model.  This  routine 

controls  all  the  forecast  time  steps  and  outputs  all  the  forecast 
fields  to  mass  storage. 

STARTF  -  Initializes  several  of  the  constants  used  by  the  model. 

FCST  -  Manages  a  single  time  step  for  the  model  forecast  variables 

(i.e.,  wind,  potential  temperature  and  pressure  of  each  sigma 
surface) .  To  take  advantage  of  the  two  central  processors  on  the 
AFGWC  UNIVAC  1110s,  the  6L  AWSPE  model  utilises  two  separate 
activities  to  perform  each  forecast.  Each  activity  forecasts 
half  of  the  model  domain  for  each  variable. 
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PACKXY  -  Packs  or  unpacks  an  input/output  (I/O)  buffer  f roa  or  into  tha 
aodel  forecast  paraaeter  arrays. 

TBNDHl  -  This  routine  coaputes  tha  tendency  of  the  potential  tesperature 
(  •  ) ,  the  pressure  thickness  (p  ),  and  the  wind  coaponents 
(  u  ,  v  )  for  the  left  half  of  the  forecast  doaain. 

TBNDH2  -  Saae  as  TBNDBl  except  it  acts  on  the  right  half  of  the  forecast 
doaain. 

TEN2H1  -  Coaputes  the  values  at  the  next  tiae  step  for  •  and  p_  for  the 
left  half  of  the  doswin.  It  also  coaputes  the  oaega  field. 

TBN2H2  -  Saae  as  TEN2H1  except  TEM2H2  acts  on  the  right  half  of  the 
forecast  doaain. 

DELPHI  -  coaputes  the  pressure-gradient-force  terns  used  to  coapute 

wind-vector  tendencies  for  the  left  half  of  the  aodel  doaain. 

DELPH2  -  Saae  as  DELPHI  except  DELPH2  acts  on  the  right  half  of  the  aodel 
doaain. 

UVMAKE  -  Coaputes  new  horisontal  wind  coaponents  froa  the  forecast  wind 
tendencies. 

TIMFIL  -  Perforas  a  tiae  saoothing  of  the  forecast  fields  for  the  wind 
coaponents,  the  pressure  thickness  fields,  and  the  potential 
teaperature. 

AVGDPl  -  Averages  the  pressure-gradient  teras  in  tiae  for  the  left  half  of 
the  aodel  doaain. 

AVGDP2  -  Saae  as  AVGDPl  except  AVGDP2  acts  on  the  right  half  of  the  aodel 
doaain. 

UVTNH1  -  Adds  the  tiae-averaged  pressure-gradient-force  contributions  to 
the  horisontal  wind  ccaponent  tendency  equations  for  the  left 
half  of  the  aodel  doaain. 

UVTHN2  -  Saae  as  UVTNHl  except  UVTHN2  acts  on  the  right  half  of  the  aodel 
doaain. 

The  above  routines  apply  the  finite-difference  equations  given  in  Section 
3  to  produce  the  desired  forecast.  Thus,  Fig.  8.1  really  represents  the 
high-level  design  of  the  6L  AtfSPB  aodel.  The  6L  AMSPB  Maintenance  Manual 
contains  a  aore  detailed  explanation. 

9.  SUMMARY 

Me  have  presented  the  aatheaatical  and  nuaerical  details  of  the  6L 
AMSPE  aodel.  Me  discussed  the  predecessor  of  the  6L  AMSPB,  the  NMC  1966  6L 


SbunuHtovtrMli  PI  aodal.  This  model  represented  •  major  stsp  forward  in 
forecast  accuracy  whan  isplemented  at  NMC  in  1966  as  wall  as  whan  implemented 
at  AIONC  in  1975. 

Lika  its  predecessor  at  MC,  the  CL  ANBPB  nodal  uses  the  hydrostatic, 
sigma-coordinate,  primitive  equations  on  a  polar  stereographic  projection  of 
the  earth.  The  finite-difference  forn  of  the  equations  is  based  on  2nd-order 
finite  difference  operators.  The  nodal  domain  consists  of  a  rectangular  area 
centered  on  the  North  Pole.  At  the  boundary  (near  ION),  the  nodal  uses  an 
isothermal,  no-slip  boundary  condition.  Terrain  is  included  in  the  model  by 
virtue  of  the  sigma  coordinate  system.  Time  and  space  smoothing  are  used  to 
control  apace  and  temporal  noise.  Ocean  boundary-layer  heat  flux  is 
parameter ised.  Surface  drag  is  applied  in  the  lowest  layer,  and 
pressure-gradient  averaging  is  used  to  allow  a  longer  tine  step. 

To  initialise  the  model,  analysed  heights  are  interpolated  linearly  in 
Bxner's  function  from  a  pressure-  to  a  sigma-coordinate  systen.  The  potential 
tamperatures  at  the  sigma  levels  are  computed.  Kinds  are  calculated  using  a 
linear  balance  relationship.  The  model  forecasts  on  the  sigma  domain  starting 
from  the  initial  conditions  using  a  centered-tine-dif ferencing  scheme.  The 
forecast  length  is  normally  72  h  for  00  and  12  GMT  data  times  and  36  h  for  the 
06  and  18  GMT  data  times.  The  forecast  fields  on  sigma  surfaces  are 
interpolated  and  stored  in  the  AfGNC  pressure-coordinate  data  base.  The 
forecast  fields  of  temperature,  wind,  height,  and  other  parameters  are 
available  for  use  by  other  computer  programs. 
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10.  INTRODUCTION 

Before  discussing  the  7-Layer  ANSPE  model,  we  will  discuss  the 
problems  that  lad  to  the  model's  development.  Me  will  than  present  some  of 
the  advantages  of  the  modal  as  compared  to  the  6-layer  version. 

Since  its  iaplementation  in  1975  the  6L  ANSPE  program  has  run  on 
AFGMC's  System  V.  However,  System  V  is  dedicated  to  satellite  data 
processing.  By  late  1979,  System  V  had  become  overloaded  with  satellite 
data.  At  that  time,  we  proposed  moving  the  6L  ANSPE  to  the  meteorological 
models  and  applications  computer.  System  I,  to  help  alleviate  this  saturation 
problem. 


afgnc  personnel  believed  after  a  new  ONIVAC  1100/81  computer  was 
installed  as  System  I,  it  would  be  possible  to  move  the  6L  ANSPE  to  System  I. 
After  several  tests,  we  demonstrated  that  the  vast  amount  of  I/O  time  required 
by  the  6L  ANSPE  made  it  impossible  for  the  model  to  meet  required  timelines  on 
the  new  System  I.  This  I/O  time  is  required  because  the  model  works  on  only  3 
vertical  strips  at  a  time.  System  I  has  a  more  powerful  CPU  than  System  V, 
but  the  I/O  processing  unit  is  not  as  efficient.  Thus,  AFGNC  decided  to 
develop  an  entirely  core-contained  version  of  the  Shuman-Hovermale  PE  model 
that  required  few  I/O  operations,  thereby  permitting  a  faster  runtime  on 
System  I. 


Improvements  made  to  the  Shuman-Hovermale  model  at  NMC  during  the 
1970's  led  T8IM  to  develop  a  core-contained  version  of  the  NMC  7L  PE  rather 
than  a  core-contained  version  of  the  6L  ANSPE.  As  stated  in  Part  I,  the  6L 
ANSPE  is  a  dry  version  of  the  original  1966  Shuman-Hovermale  6L  PE  model. 

Some  modifications  were  made  during  the  mid  1970's  to  the  6L  ANSPE  to  allow  a 
longer  time  stepi  pressure-gradient  averaging,  tendency  truncation  south  of 
20N,  and  horisontal  space  smoothing.  These  changes  did  not  significantly 
improve  forecast  accuracy.  Thus,  AFGNC  forecasters  and  others  (e.g.,  Leary, 
1971)  continued  to  identify  problem  areas  in  the  6L  ANSPE  forecasts.  Two 
examples  are  the  slow  phase  speed  of  the  shorter,  weather  producing  waves  and 
the  generation  of  spurious  vorticity  patterns  near  the  model  boundaries.  NMC 
solved  both  of  these  problem  areas  in  the  7L  PE  model  by  incorporating  better 
physical  and  numerical  techniques. 


To  reaedy  som  of  tb«  CL  AMSPB  aodel  disadvantages  and  to  efficiently 
take  advantage  of  the  large  aeaory  available  (1,000,000  words)  on  System  I, 
Krone  decided  to  develop  2nd-  and  4th-order  core  contained  versions  of  the  me 
7L  PE  developed  by  Caapana  (1977,  1978,  1979a,  and  1979b).  APGNC/DO  tasked 
TSIM  with  this  technology  transfer  project.  The  large  aaount  of  aeaory  this 
aodel  uses  (400,000  words)  eliainates  90  percent  of  all  I/O  operations  and 
should  allow  the  7L  AN8P1  to  run  operationally  on  Systea  I.  Furtheraore,  the 
iaproved  nuaerical  forecast  techniques  (4th  order  differencing  and  better 
boundary  conditions)  will  definitely  iaprove  the  aodel  forecasts. 

10*1  The  7-Layer  AMSPB  Model. 

When  coaparing  the  7L  and  6L  ANSPE  aodels,  the  aajor  iaproveaents  in 
the  7L  version  lie  in  the  aore  sophisticated  nuaerical  techniques  used.  The 
specific  techniques  are  listed  below. 

(1)  The  7L  ANSPE  has  the  capability  to  use  4th-order  finite 
differencing  for  the  advection  teraa  in  the  equations  of  notion.  Caapana 
(1978)  showed  that  the  1-bedient  grid  version  of  the  NNC  71  PE  aodel  with 
4th-order  finite  differencing  verifies  nearly  as  accurately  as  the  2nd-order, 
1/2-bed lent  version.  Both  the  4th-order,  1-bedient  and  2nd-order,  1/2-bedient 
aodels  produce  better  forecasts  because  they  handle  the  propagation  speeds  of 
the  shorter  weather  producing  waves  aore  accurately  than  a  1-bedient, 

2nd-order  aodel.  The  coaputer  cost  of  a  1/2-bedient  aodel  is  approximately  8 
tines  greater  than  for  a  1-bedient  aodel,  while  only  about  a  40  percent  cost 
increase  is  involved  in  going  froa  2nd-  to  4th-order  finite  differencing  at 
the  sane  grid  increaent.  Caapana  (1978)  showed  that  4th-order  differencing 
was  a  very  cost-effective  way  to  iaprove  nuaerical  forecasts. 

(2)  The  7L  ANSPE  uses  the  "energy-vorticity"  fora  for  the 
equations  of  notion  while  the  CL  ANSPE  uses  the  "standard”  u  and  v  wind 
coaponent  fora  of  the  equations  of  notion.  Shuaan  and  Stackpole  (1968)  showed 
that  the  finite-difference  fora  of  the  "standard”  set  of  equations  leads  to 
spurious  vorticity  production,  while  the  "energy-vorticity”  finite-difference 
fora  does  not.  Forecasters  also  have  noticed  a  significant  aaount  of  spurious 
vorticity  generated  at  the  boundaries  in  the  current  6L  ANSPE  aodel. 

(3)  The  7L  ANSPE  forecasts  winds  and  potential  teaperatures  at  7 
layers  in  the  aodel  ataosphere.  The  6L  ANSPE  forecasts  the  wind  and 
potential  temperature  for  6  layers.  MC  has  found  that  by  asking  the 
additional  layer  an  active  forecast  layer,  aodel  verification  improves 
slightly  in  the  upper  layers. 

(4)  The  7L  ANSPE  uses  energy-conserving  finite-difference  schemes 
for  computing  geopotential  heights  froa  the  hydrostatic  equation  and  for 
applying  dry  convective  adjustment.  Arakawa's  (1972)  energy-conserving 
scheme ,  as  given  by  Phillips  (1974),  is  used  to  compute  the  geopotential 
height. 

(5)  Other  differences  between  the  CL  and  7L  ANSPE  aodels  do 
exist.  These  will  be  covered  by  the  aore  detailed  description  of  the  7L  ANSPE 
aodel  given  in  the  following  sections.  Section  11  will  cover  the  model 


primitive  equations,  while  Section  12  will  present  the  finite-difference 
formulation  of  the  7L  AN8PB.  The  marching  process  will  be  outlined  in  Section 
13.  Sections  14  and  15  will  discuss  the  initialisation  and  output  fields, 
respectively.  In  Section  16  we  will  present  the  7L  AWSPE  model  flow  diagram 
and  abridged  maintenance  manual. 

11.  THE  7-LAYBR  ASMPE  MODEL  CHARACTERISTICS 

11.1  The  Primitive  Equations. 

Caapana  (1978)  presented  the  energy-vorticity  form  of  the  hydrostatic 
meteorological  equations  for  a  model  with  a  vertical  sigma  coordinate  and 
horisontal  Cartesian  coordinates  on  a  conformal  polar  stereographic  projection 
of  the  earth  as: 
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where  E  is  defined  ast 
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(11.6) 


E  *  %(uj  *  v*)* 

where 

u0  *  uw»  v0  *  Vm  *  (11.7) 


These  six  equations  are  written  in  the  six  dependent  variables  ur  v,  8,tt ,  <J>, 
and  pg.  The  variables,  u  and  v,  are  the  horisontal  wind  components,  e  is  the 
potential  temperature,  f  is  the  Coriolis  parameter,  a  is  the  aap  factor,  9  is 
the  geopotential  height,  p0  is  the  pressure  thickness  between  the  three  sigma 
domains)  and  it  is  Exner's  function.  Q  represents  heating  and  for  the  current 
version  Q  represents  only  the  sea  surface  heat  flux.  F  represents  friction 
which  is  only  used  in  the  7L  MtSPE  in  the  fora  of  surface  drag.  The 
overriding  bar,  (  ) ,  in  (11. 5)  represents  a  vertical  average  within  each 

sigma  domain.  The  (j  in  (11.1)  and  (11.2)  represents  the  pressure  gradient 
averaging  technique  given  by  Brown  and  Caapana  (1978) . 

11.2  The  Vertical  Domain. 

Fig.  11.1  is  a  depiction  of  the  vertical  structure  of  the  7L  AMSPE 
model.  The  model  contains  seven  vertical  layers  which  are  partitioned  into 
three  separate  sigaa  domains.  The  lowest  layer  is  a  planetary  boundary  layer 
with  the  vertical  coordinate  defined  by 


The  second  sigma  domain  consists  of  three  layers  and  covers  the  remainder 
of  the  troposphere.  This  domain  is  defined  by 

p-p** 

°T  ’  |p,-W|.p~  «“•»» 

where  p**  is  the  pressure  at  the  tropopause.  The  tropopauae  is  considered  to 
be  a  material  surface)  i.e.,  ®3  ■  0. 

The  third  sigma  domain  consists  of  three  layers  and  constitutes  the  model 
stratosphere.  The  stratospheric  domain  is  defined  by 
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(11.10) 


The  top  o£  the  nodel  is  specified  at  p  ■  50  mb. 


11.3  The  Horizontal  Grid. 


The  7 L  AHSPE  uses  the  horizontal  grid  described  in  Section  2.3. 

11.4  Prognostic  and  Diagnostic  Variables. 

The  forecast  or  history  variables  are  u,  v,  pg,  and  e  .  The  diagnosed 
variables  are  it  ,  z,  o  ,  and  p.  In  Pig.  11.1,  we  show  that  u,  v,  z,  p  and  -n 
are  defined  at  the  half  levels  and  a  and  p  are  defined  at  the  full  levels. 
Actually,  a  is  defined  at  the  center  point  of  each  square  formed  by  four 
full-level  grid  points  (see  Fig.  11.2).  The  methods  used  to  calculate  these 
variables  will  be  presented  in  the  next  two  sections. 

12.  SECOND-  AND  FOURTH-ORDER  FINITE-DIFFERENCE  FORMULATION  OF  THE  7-LAYER 
AWSPE  MODEL. 

In  this  section,  we  present  the  second  and  fourth-order  finite- 
difference  formulations  of  the  model,  the  parameterization  of  planetary 
boundary  layer  (PBL)  friction,  the  method  of  calculation  of  6  ,  the 

lateral  boundary  conditions,  the  horizontal  smoothing,  the  time  smoothing,  the 
pressure-gradient-averaging  methods,  and  the  technique  used  to  compute  ui  . 

12.1  The  7-Layer  Second-  and  Fourth-Order  AWSPE  Finite-Difference 
Equations. 

Again,  we  will  use  the  shorthand  notation  given  in  Part  I  (Shuman  and 
Hovermale,  1968).  The  operator  used  to  construct  the  second-order  version 
(i.e.,  P{  and  fx)  will  remain  the  same.  Campana  (1977  introduced  the 
fourth-order  operators  using  similar  notation.  He  defines  fX|,  and  P*hto 
be  fourth-order  operators  for  the  derivative  and  the  sum  as: 
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As  in  ths  sscond-ordsr  formulation,  the  derivative  and  averaging 
operator  are  valid  Midway  between  a  set  of  four  grid  points.  Thus,  the 
derivative  is  used  to  go  to  the  half  grid  points  and  the  sub  is  used  to  return 
to  the  whole  grid  point. 

Now  define  the  fourth-order  operators  for  the  horisontal.  First,  it 
is  instructive  to  define  the  second-order  operators  for  coaparison  purposes. 
These  operators  for  a  horisontal  (x,y)  coordinate  system  are  given  byt 


**  ”  2Ax  +  <12‘3) 
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***  “  fi*Hj,j-%  +  -h  *  fiJj,W 
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where  i  refers  to  grid  points  in  the  x-direction  and  j  refers  to  grid  points 
in  the  y-direction.  The  fourth-order  operators  are  found  by  combining  (12.1) 
and  (12.2).  For  a  horisontal  (x,y)  coordinate  system,  these  operators  are 
given  by 
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(12.7) 
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9^f±+%,j+if  +  fi+34,j-%  +  fi-%,j-Hj  +  fi-34,3J* 

+  fi-Hj,j+%+  fi-Hl,j-34  +  fl-Js.j+%  +  fl-J5,j-3/2) 
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These  fourth-order  operators  require  nany  aore  calculations  than  the 
cosqparable  second-order  ones  given  by  (12.3)  through  (12.5). 

Next,  we  present  the  second-  and  fourth-order  versions  of  the  7L  AWSPE. 
Caepana  (1978)  found  that  equally  accurate  forecasts  could  be  obtained  either 
by  (1)  using  the  fourth-order  operators  for  the  non  linear  advection  teres  and 
final  overbar  operation  only  or  by  (2)  using  the  fourth-order  operators  for 
every  operation.  Therefore*  the  fourth-order  version  we  developed  uses 
fourth-order  differencing  only  for  the  nonlinear  advection  teres  and  the  final 
overbar  operation. 

Using  the  Shuean-Hovereale  finite-difference  notation  and  neglecting 
heating  and  friction*  the  finite-difference  representations  of  (11.1)  through 
(11.5)  for  the  second-order  finite-difference  nodel  can  be  written  (Caapana* 
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where 
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The  l  I  over  the  pressure-gradient-foree  terse  in  (12.9)  and  (12.10) 
represents  the  pressure-gradient-averaging  technique  of  Brown  and  Caapana 
(1978). 

The  corresponding  fourth-order  equations  are  written  as: 


|S  -  ♦  T*  |  S»*  -  i>li  -  t»k  -  &  ♦  cWih  -  li  5 *»)' 

3*  xh  yh  xh  x  p  X  ° 

|£  »  -  u^tf^  +  -  axW]  -  "  (oV^J 


(12.13) 


12.14) 


II  .  V  fiWfc  -  (66*^) 

3-c  a 


a1  ■  *  <‘oi‘'sx%)]  '  £lV4 

X«  Ui.  * 


(12.15) 


(12.16) 


where  U  ■  u;  V  ■  v; 


fa1  '* 


Note,  as  stated  before,  fourth-order  differencing  is  used  only  for  the 
nonlinear  advection  teres  and  for  the  final  averaging  (final  overbar) 
operation. 

12.2  PBL  Friction. 

For  the  7L  AMSPE,  frictional  dissipation  is  considered  only  in  the 
model  boundary  layer.  The  formulation  of  this  frictional  dissipation  is  the 
same  as  given  in  Part  I,  Section  3.2. 

12.3  The  Calculation  of  "Sioma  Pot*. 

In  a  manner  similar  to  the  6L  AMSPE  model,  the  method  used  to  compute  o 
for  the  7L  AMSPB  follows  almost  directly  from  the  finite-difference  fora  of 
the  pg  tendency  equation  (3.11)  and  the  fact  that  o  ■  0  at  levels  a  ■  7*, 
3**,  and  0. 

First,  po  ■  50  mb  in  the  boundary  layer  and  is  constant  in  tine  and 
space.  Thus,  if  we  use  the  non-vertically  integrated  fora  of  (12.16)  and 
remove  the  (  )*/  we  have 
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where 
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and 
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A  different  approach  suit  be  employed  in  the  troposphere  since  P0  varies 
both  in  tine  and  space.  However,  Pa  is  linear  with  respect  to  a,. 

Therefore  Paa  is  equal  to  sero.  Using  this  fact,  an  equation  for  °aa  can 
be  found  and  then  (0J5  and  (6^)4  can  be  computed.  Differentiating  (3.11) 
with  respect  to  a  yields 
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Thus  defining 
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(12.21)  can  be  written  as 
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Applying  (12.23)  at  levels  k  ■  5  and  4  yields 
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Solving  (12.24)  foe  (oT)4  and  (aT) 5  and  given  that 
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and  (ot)3  *  0  yields 
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Similarly,  in  the  stratosphere 
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12.4  Lateral  Boundary  Condition*. 

The  7L  AWSPE  has  the  capability  to  use  either  the  no-slip  boundary 
conditions  described  in  Section  3.4,  or  the  sponge  boundary  conditions 
developed  by  Deaven  at  MMC  in  1979.  We  will  outline  the  sponge  boundary 
conditions  (me  Numerical  Weather  Prediction  Activities  Report,  1979) . 

The  sponge  boundary  is  designed  to  provide  a  better  forecast  near  the 
outer  edges  of  the  model  dosutin  by  dasping  waves  traveling  toward  or  away  from 
the  boundary.  To  daap  these  waves,  the  following  term  is  added  to  the  u,  v, 

6,  and  Pa  tendency  equations  at  the  second  through  sixth  grid  rows  fro*  the 
boundary! 
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where  ♦  represents  u,  v,  P0  and  •  i  4°  is  the  value  at  the  initial  time; 
and  K  is  0.04  for  the  7 L  AWSPE.  The  vslue  at  the  boundary  (grid  row  1)  is 
maintained  at  the  initial  value  throughout  the  forecast  period.  Thus,  the 
rows  near  the  boundary  are  continually  nudged  towards  the  initial  boundary 
values. 

12.5  Thermodynamic  Forcing 

Thermodynamic  forcing  in  7L  AWSPE  model  is  the  same  as  the  (L  AWSPE 
(Part  1,  Section  3.5).  As  before,  the  only  forcing  is  due  to  the  eddy  flux  of 
heat  frosi  the  ocean  surface. 


12.6  Horliontal  Smooth! 


The  7L  AWSPE  model  can  use  two  additional  forms  of  smoothing  besides 
the  snoothing  included  in  the  finite-difference  formulation  of  the  model.  The 
first#  which  is  always  used#  is  the  diffusive  space  smoother  given  by  Shuaen 
(1977) .  The  plan  of  this  smoother  is 
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where  |i  is  chosen  as  0.98  for  the  7L  AWSPB.  That  is#  each  variable  is 
smoothed  by  applying  the  above  plan  at  each  grid  point.  The  weight  at  the 
central  grid  point  where  the  smoothing  result  is  added  to  the  original  value 
is  given  by  y2_1  .  Other  weights  are  as  specified  in  the  plan. 

The  second  fora  is  applied  below  20N.  mis  smoother  is  the  tendency 
truncation  is  given  in  Section  3.6  of  Part  1. 

12.7  Time  Smoothing. 

Because  the  model  uses  leapfrog  time  differencing  scheme,  a  tiae 
smoother  is  applied  to  control  tiae  splitting.  The  tiae  smoother  has  the  fora: 


-  2/  ♦  *n+T)  (12.32) 


where  q*  represents  the  variables  u#  v#  •#  and  pQ  ,  n  represents  a  given  tiae 
level#  and  a  ■  0.075. 

12.8  Computation  of  Omega 

Omega  (u>)  is  diagnosed  for  use  by  forecasters  and  applications  programs. 
The  model  does  not  use  the  m  value  in  any  way*  "Sigma  dot"  ( q  )  is  the 
model's  vertical  aotion  parameter. 

is  computed  directly  using 
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or  in  finite  difference  fora 
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The  value  of  u>  is  then  averaged  over  a  2-h  period  (€  tiae  steps  for 
which  At  *  1200  s)  before  output.  For  instance,  the  value  "valid"  at  the 
3-h  forecast  point  would  be  an  average  froa  1+20  to  3+20. 

12.9  Coaputatlon  of  the  Ceopotential  Height 

Unlike  the  6L  AMSPB,  the  7L  AMSPB  does  not  use  the  siaple 
finite-difference  fora  of  the  hydrostatic  equation  (i.e.  (3.9))  used  by  the  6L 
AMSPB  to  coapute  the  geopotential  height.  Instead  the  7L  AMSPB  uses  a 
technique  very  siailar  to  the  one  given  by  Phillips  (1974) .  Phillips 
technique  was  aodeled  after  Arakawa's  (1972)  energy  conserving  vertical 
integration  technique.  Me  will  only  suaaarise  the  basic  technique  here. 
Phillips  (1974)  and  Brown  (1974)  gave  a  aore  detailed  description.  First  we 
define 
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where  A^  is  the  separation  between  sigaa  levels  k  and  k+1,  see  Fig.  (12.1). 
Nest  define 
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where  the  hat  (A  )  denotes  a  variable  at  the  whole  level,  while  the  absence  of 
the  (*)  represents  a  variable  at  the  interaediate  levels  and puis  pressure  at 
the  top  of  a  given  signa  doaain.  Mow  the  geopotential  can  be  ccaputed  as 
follows  t 
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Figure  12.1.  Generalized  depiction  of  variable  location  for  computing  the 
geopotential  height  using  the  method  given  by  Phillips  (1974) 
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b.  Compute  Bxner's  function  for  the  whole  levels. 
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c.  Compute  ir^  fro* 
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d.  Coopute  predicted  values  of  6  froa  the  predicted  values  of  6  and  e 
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where  is  the  surface  value  of  the  geopotential. 

f.  Finally,  determine  the  geopotential,  <f>t  ,  k-2,3,...,K  from 
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The  above  technique  is  used  in  7 L  AWSPE  routine  GEOPTL  which  was  adapted 
from  NMC's  7L  PE  routine  GE0P7L. 

12.10  Dry  Convective  Adjustment. 

The  7L  AWSPE  model  applies  a  dry  convective  adjustment  procedure  to 
maintain  model  stability.  If  this  procedure  is  not  used,  the  phase  speed  of 
internal  gravity  waves  becomes  unstable  (complex)  when  superadiabitic 
conditions  exist. 

The  procedure  used  to  apply  the  dry  convective  adjustment  was  taken 
directly  from  the  NMC  7L  PE  model  software.  A  similar  procedure  is  outlined 
in  Haltiner  and  Williams  (1980)  and  will  be  repeated  here.  Both  of  the 
schemes  conserve  total  potential  energy;  that  is. 


STpdz  ■  /  fPt  JTdp  •  0  (12-43> 
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where  T  and  B  represent  the  top  and  bottom  of  the  unstable  layer, 
respectively.  Physically,  dry  convection  develops  whenever  the  lapse  rate 
exceeds  the  dry  adiabatic  value.  This  transports  heat  upward  until  a  neutral 
lapse  rate  results.  Thus,  potential  energy  is  converted  to  kinetic  energy 
which  is  eventually  converted  into  heat.  The  result  is  a  total  energy 
redistribution  of  potential  temperature. 

The  71.  AWSPE  essentially  uses  a  formulation  similar  to  (12.43)  to  solve 
for  the  energy  conserving  readjustment  of  temperature.  For  one  unstable  layer 
between  levels  k  and  k+1,  (12.43)  can  be  written  as 
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Thus  (12.44)  and  (12.45)  can  be  solved  for  and  to  gefc  the  dry 

convective  adjustment  to  the  temperature  field.  A  very  similar  technique  is 
applied  at  all  levels  in  the  7L  AWSPE  model  to  ensure  no  superadiabatic  lapse 
rates  are  allowed  to  exist. 

12.11  Fourth-Order  Versus  Second-Order  Differencing. 

The  7L  AWSPE  can  use  either  fourth-order  or  second-order  finite 
differencing  for  the  advection  terms  in  the  equations  of  motion.  Campana 
(1978)  found  that  fourth-order  differencing  greatly  improved  forecast 
accuracy.  Table  12.1  presents  his  results  and  demonstrates  the  expected 
improvement  in  the  fourth-order  version  of  the  71.  AWSPE  as  compared  to  the 
second-order  version. 

Table  12.1  contains  the  500-mb  RMS  vector  wind  error  and  500-mb  mean 
RMS  temperature  error  at  48h.  Note  that  the  fourth-order  1-bedient  7L  PE 
verifies  almost  as  well  as  the  (then  operational)  second-order  half-bedient 
NMC  7L  PE.  The  improved  forecast  accuracy  of  the  fourth-order  NMC  7L  PE  (a 
moist  version  of  the  7L  AWSPE)  is  impressive  when  compared  to  the  NMC  6L  PE  (a 
moist  version  of  the  6L  AWSPE) . 

Haltiner  and  williams  (1980)  showed  fourth-order  differencing  improves 
phase  speed  accuracy  for  waves  whose  length  is  4  to  12  grid  increments  (DX) 
for  the  simple  advection  equation.  Here  DX  refers  to  the  separation  between 
adjacent  grid  points  in  the  model  (1-bedient  for  the  6L  and  7L  AWSPE  models) . 
Thus,  the  improved  forecast  accuracy  found  in  Campana 's  (1978)  fourth-order 
results  was  due  to  improved  forecasts  of  phase  speeds  for  wavelengths  of  four 
grid  increments  and  larger. 


Table  12.1.  Fourth-order  versus  second-order  finite  differencing  forecast 
accuracies  for  selected  cases  (Caapana,  1978) 


MEAN  500-ab 
RMS  Vector 
Wind  Error 

(48-h) 

(■  s"1) 

NMC 

6L  PE 

(1-bed ient) 

OPNL 

Second  Order  NMC 

7L  PE 

(1/2-bedient) 

Fourth  Order  NMC 
7L  PE 

(1-bedient) 

1000  ab 

9.53 

7.67 

7.82 

500  wib 

11.13 

8.66 

8.89 

300  mb 

16.26 

12.77 

13.04 

100  mb 

8.13 

8.27 

Mean  (48-h)  RMS 

500-mb  Temperature  Error 

(°C) 

lOOOab 

3.78 

3.13 

3.16 

SOOeb 

3.06 

2.51 

2.56 

300»b 

3.08 

2.57 

2.59 

lOOnb 

— 

3.48 

3.53 

13.  THE  7-LAYER  AMSPE  MODEL  MARCHING  PROCESS 


The  finite-difference  equations  for  the  second-order  Model  are  solved 
in  the  following  manner.  For  the  fourth-order  model  replace  (12.9)  through 

(12.12)  with  (12.13)  through  (12.16).  Given  initial  values  for  u,  v,  e, 
and  Pa, future  values  of  these  variables  are  computed  using  (12.9)  through 

(12.12) .  To  apply  pressure-gradient  averaging,  the  future  values  of  PQ  and  e 
must  be  computed  before  the  pressure-gradient  force  can  be  calculated. 

Thus,  o  is  computed  using  the  method  given  in  Section  12.3.  Then  Pa  is 
forecast  from  (12.12),  and  e  is  computed  using  (12.11).  The  pressure  at  all 
levels  is  computed  for  the  future  time  step  N+l  (where  N  represents  an 
arbitrary  time  step)  by  summing  PQ  downward  from  the  top  level  (0) ,  where  p  - 
50  mb.  Now  Exner's  function,  it  ,  at  N+l  is  computed  from 
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where  pQ  ■  1000  mb.  The  technique  given  in  Section  12.9  is  then  used  to 
compute  the  heights,  z,  at  the  intermediate  levels  for  time  N+l  given  the 
terrain  height  at  level  7.  Finally,  the  pressure-gradient  force  can  be 
computed  from  the  values  of  z  and  -nr  at  times  (N-l,  N,  N+l)  using  (3.12). 
Then  (12.9)  and  (12.10)  are  solved  for  the  values  of  u  and  v  at  time  step 
N+l.  He  repeat  this  process  in  time  until  the  desired  forecast  length  is 
reached.  After  the  first  time  step,  a  centered  finite  difference  scheme  is 
used  for  the  time  derivative.  Because  pressure-gradient  averaging  is  used  in 
the  model,  a  time  step  of  1200s  is  permitted  for  the  second-order  version. 

For  the  fourth-order  version  the  time  step  must  be  reduced  to  900s.  For  the 
first  time  step,  both  models  use  a  forward-finite-difference  scheme.  A 
forward  scheme  is  also  used  for  surface  drag  and  sea  surface  heating  terms. 


14.  INITIALIZATION. 

The  initialisation  process  is  detailed  in  Shuman  and  Hovermale 
(1968).  This  process  for  the  7L  AMSPE  is  identical  to  that  given  in  Section 
5,  Part  I,  with  two  exceptions.  First,  data  from  the  APGMC  High  Analysis 
Model  (HIANAL)  at  50  mb  are  used  to  compute  values  for  the  top  layer  of  the 
model.  Values  in  the  top  layer  of  6L  AMSPE  model  are  chosen  in  a  more 
arbitrary  manner.  Second,  the  nonlinear  balance  equation  is  used  vice  the 
linear  balance  equation  in  Part  I. 

15.  OUTPUT  FIELDS. 

The  7L  AMSPE  model  produces  the  same  pressure-coordinate  data  base  as 
the  6L  AMSPE  (Section  6) .  The  technique  used  to  convert  7L  AMSPE 
sigma-coordinate  forecast  fields  to  the  pressure-coordinate  data  base  is  also 
identical  to  that  used  by  the  6L  AMSPE  (Section  6) . 


16.  A  DESCRIPTION  OF  THE  SUBROUTINES  AND  MAJOR  ARRAYS  USED  IN  THE  7-LAYER 
AWSPE  PROGRAM. 

In  this  section  we  describe  the  subroutines  used  in  the  7-Layer  AMSPE 
prograa.  Then  we  present  a  description  of  the  major  common  blocks  and  arrays 
used  throughout  the  model.  This  section  is  intended  for  those  requiring  a 
more  detailed  knowledge  of  the  7L  AWSPE.  Further  explanation  of  the  actual 
details  of  the  model  can  be  found  in  the  in-line  program  documentation.  Figs. 
16.1/  16.2/  and  16.3  give  the  7L  AWSPE  model  subroutine  flow  diagram.  The 
numbers  represent  the  order  in  which  the  individual  routines  are  called.  The 
basic  functions  of  each  subroutine  are  given  below. 

PEFCST-  The  main  routine  or  driver  for  the  7L  AWSPE  model.  It  sets 

constants/  reads  SSETC  information,  inputs  the  initial  data  on  sigma 
coordinates  from  absolute  PBINIT,  produces  the  desired  forecast,  and 
outputs  the  forecasts  to  files  for  conversion  to  the 
pressure-coordinate  data  base. 

BALNCE-  Computes  the  u  and  v  components  of  the  wind  from  the  sigma  surface 
geopotential  using  the  nonlinear  balance  equation. 

BALUV-  Computes  the  wind  from  the  height  and  Exner's  function  on  a  sigma 
surface  using  the  nonlinear  balance  equation. 

BNDRY-  Computes  the  u  and  v  boundary  conditions  by  applying  the  "diffusive 
nudge"  boundary  conditions.  (Section  12.4). 

btdry-  Computes  the  6  and  Pa  boundary  conditions  by  applying  the 
"diffusive  nudge”  boundary  conditions.  (Section  12.4.) 

BNDSET—  Computes  the  Laplacian  of  the  the  initial  boundary  values  (6  grid 
rows) .  These  values  are  used  by  subroutine  BNDRY  to  compute  the 
"diffusive  nudge”  lateral  boundary  conditions.  (Section  12.4). 

DIAG  -  Computes  all  diagnostic  energetics. 

Div-  Computes  the  divergence  on  the  earth.  This  diagnostic  information  is 
not  used  by  the  6L  AWSPE  model. 

DRY-  Applies  dry  convective  adjustment  to  the  forecast  potential 

temperature  field,  if  6  above  is  less  than  e  below,  then  the 
scheme  adjusts  the  two  potential  temperatures  to  produce  a  neutral 
lapse  rate.  (See  Section  12.10) 

FCST-  Forecast  module  driver.  This  routine  forecasts  future  values  of  u,  e 
v,  and  Pa  using  centered  time  differencing  with  a  time  step  of  2 
*DT  seconds  (where  DT-1200  s  for  second  order  and  DT  ■  900  s  for 
fourth  order) .  Forecast  is  computed  on  a  polar  stereographic 
projection. 


Figure  16.1.  The  subroutine  flow  diagram  for  the 
Lower  Portion  of  the  7L  AHSPE  Model . 
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PCSTUV-  Given  the  tendency  values  computed  by  UVPCST,  it  forecasts  u  and  v  at 
the  next  tiae  step.  It  then  adjusts  the  boundary  values  and  applies 
space  and  tiae  saoothing. 

PILERR-  This  is  an  1/0  error  teraination  routine. 

PIXARR-  Switches  the  overlap  areas  of  each  Data-Bank  for  the  potential 

teaperature  array  (THETA) .  That  is,  it  stores  rows  27  and  28  from 
D-Bank  1  into  rows  31  and  32  in  D-Bank  2.  It  then  store  rows  33  and 
34  from  D-Bank  2  into  rows  29  and  30  in  D-Bank  1. 

PIXCON-  Switches  the  overlap  areas  of  each  D-Bank  for  the  land-sea  indicator, 
the  sea  surface  teaperature,  the  drag  coefficient,  the  map  factor, 
and  the  Coriolis  parameter.  Sane  as  PIXARR  for  the  sea  surface 
temperatures,  the  drag  coefficient,  the  nap  factor,  and  the  Coriolis 
parameters. 

PIXPHI-  Same  as  PIXARR  only  for  the  geopotential  height  field. 

PIXWRK-  Transfers  values  from  the  left  side  (J-33,  34)  of  that  portion  of  the 
work  array  where  J*31,  61  to  the  right  side  (J*29,  30)  of  the  work 
array  where  >1,  30.  Values  are  then  transferred  from  the  right  side 
(J*27,28)  of  that  portion  of  that  work  array  where  J«l, 30  to  the  left 
side  (J>31,  32)  of  that  portion  of  the  work  array  where  J»30,61. 

PIXUV-  Switches  the  overlap  areas  of  each  Data-Bank  for  the  u  and  v  wind 

component  arrays.  It  stores  rows  27  and  29  from  D-Bank  1  into  rows 
31  and  32  in  D-Bank  2.  It  then  stores  rows  33  and  34  from  D-Bank  2 
into  rows  29  and  30  in  D-Bank  1. 

FXPSIG-  Same  as  PIXARR,  except  for  the  array. 

FXSDOT-  Same  as  PIXARR,  except  for  the  'a  array. 

GEOPTL-  Computes  the  geopotential  height  field  at  all  levels  for  a  given 

time.  The  method  used  here  follows  from  Brown  (1974)  and  Phillips 
(1974),  see  Part  II,  Section  12.9.  This  method  was  first  developed 
by  Arakawa  and  an  example  of  the  procedure  is  given  in  Arakawa 
(1977).  The  code  used  here  follows  directly  from  NMC's  subroutine 
GB0P7L. 

GBOPTL1-  Same  as  GEOPTL  except  it  computes  the  geopotential  height  for  a 
single  grid  point. 

GB0PTL2-  Computes  the  geopotential  height  field  at  all  levels  for  a  given 
time.  (Section  3.9). 

GRAPH-  Diagnostic  program  that  prepares  a  graph  of  the  7L  AHSPE  model 
kinetic  energy  amplification  with  respect  to  time.  This  is  a 
diagnostic  routine  only. 


HEAT- 


INITIAL- 


INITL2- 

KEPCAL- 

LDBNK1- 

LDBNK2- 

NFM02- 

OBDSET- 

OLBDTH- 

OLDBRY- 

OMEGA- 


Coaputes  the  only  diabetic  tera  used  by  the  7L  awspe.  Sea  surface 
heating  is  parameterized  using  the  specification  given  by  Shuman  and 
Hovermale  (1968).  Other  diabetic  forcing  teras  could  be  added  here. 
(Section  12.5). 

Reads  the  initial  data  file  produced  by  absolute  PEINIT.  For  a 
noraal  6XQT  the  initial  data  is  stored  on  file  10.  Otherwise  files 
L3HR  and  99  are  used  as  the  current  and  current  minus  one  time  steps 
in  order  to  start  with  leapfrog  tiae  differencing.  The  initial 
fields  are  on  the  sigaa  coordinate  systea  shown  in  subroutine  PCST. 
The  routine  does  the  following:  First,  it  deteraines  whether  this  is 
a  restart  or  a  forecast  starting  from  tiae  0.  Next,  the  data  are 
retrieved  and  stored  in  the  proper  arrays.  The  subroutine  calls 
GBOPTL  to  confute  the  height  field  in  preparation  for  the  first  time 
step.  Finally,  diagnostic  output  fields  are  generated,  if  requested. 

Initializes  the  model  with  Rossby-Haurwitz  waves  to  test  the  model 
numerics.  Many  other  initialisations  are  available  in  this  test 
routine.  This  is  only  a  test  routine. 

Computes  the  product  of  the  kinetic  energy  and  the  aap  factor  for  use 
in  UVFCST.  Place  the  result  in  B  (I,  1,  1,  k) .  Computes  the 
pressure  gradient  force  for  UVFCST. 

Sets  the  D-Bank  pointer  or  Bank  Descriptor  Index  (BDI)  to  D-Bank  1. 
The  arrays  in  coaaon  blocks  FIKLD1  and  FIBLD2  are  mirror  iaages  of 
one  another.  Common  block  FIBLD1  contains  all  arrays  that  are  mapped 
into  D-Bank  1  and  cowon  block  FIKLD2  contains  all  arrays  that  are 
mapped  into  D-Bank  2.  The  array  names  in  FIELD2  are  treated  as  dummy 
naaes  and  are  never  referenced  by  the  model,  but  the  storage 
locations  are  used  by  changing  the  Bank  Descriptor  Indicies  (BDI) 
through  calls  to  LDBMR1  and  LDBNR2. 

Same  as  LDBNK1  except  it  sets  the  BDI  or  D-Bank  pointer  to  D-Bank  2. 

Smooths  the  wind  field  after  it  is  nonlinearly  balanced  to  the  height 
field. 

Applies  the  old  6L  AMSPB  model  lateral  boundary  conditions  to  e  and  p 
(Section  3.4).  This  routine  is  used  during  the  initialisation  0 

process. 

Saae  as  0BD8BT  except  this  routine  is  used  during  the  7L  AMSPB  model 
marching  process. 

Applies  the  old  6L  AMSPB  model  lateral  boundary  conditions.  These 
conditions  are  applied  to  u  and  v  in  this  routine.  (Section  3.4). 

Computes  OMBGA  for  use  by  applications  programs  and  AF6WC/MF 
forecasters.  This  is  a  diagnostic  routine  only.  See  section  12.8. 


OUT-  Display*  a  53  X  57  array.  One  half  of  the  array  will  be  displayed 
per  page.  If  less  than  7  levels  are  requested,  the  routine  outputs 
starting  from  the  bottoai  level. 

OUTPUT-  Checks  to  see  if  it  is  time  to  store  into  the  data  base  files.  Every 
three  hours,  data  will  be  stored  in  the  files.  Every  12  hours, 
either  PE12,  PE24,  PE36,  or  PE48  will  be  started  to  store  the  data 
into  the  pressure-coordinate  data  base. 

PICAL-  Computes  Exner's  function  at  the  intermediate  vertical  levels  using 
(12.37).  This  technique  helps  to  conserve  energy  in  the  calculation 
of  the  pressure  gradient  force. 

PRES-  Computes  the  pressure  of  each  point  on  the  sigma  surface. 

PSGFCT-  Forecasts  Pc  at  the  next  time  step  from  the  tendency  value  which  is 
computed  by  PSIGMA.  The  boundary  values  are  adjusted  and  the 
forecast  values  are  time  smoothed. 

PSIGMA-  Computes  the  tendency  of  Pa  for  the  model  troposphere  and 

stratosphere.  Pa  is  a  constant  (SOmb)  in  the  (one-layer)  boundary 
layer.  To  compute  pQ  ,  we  use  (12.12)  or  (12.16)  from  Section  12.1. 

SHTBDY-  Applies  the  old  6L  AWSPE  model  boundary  conditions.  Here,  these 

conditions  are  only  applied  to  the  terrain  heights  to  prevent  noise 
on  the  boundaries  caused  by  variable  terrain  at  the  boundary.  See 
Section  3.4. 

SIGDOT-  computes  the  sigma  coordinate  vertical  motion,  o  at  time  "N”  for  use 
in  the  tendency  equations.  (Section  12.3.) 

SMOOTH-  Smooths  array  ANK1  (time  level  H  -1)  and  then  adds  the  result  to 
array  AMP1  (time  level  H+l) .  Space  smooths  the  6  and  p  fields. 
(Section  12.6).  0 

SMTHUV-  Same  as  SMOOTH  except  for  u  and  v. 

START-  This  routine  is  called  by  main  program  PEFCST  to  retrieve  the  0SETC 
information.  The  many  ESBTC  options  available  are  given  in  Table 
16.1. 

STARTF-  Sets  the  basic  parameters  for  the  7L  AWSPE  model,  e.g.,  time  step, 
gas  constants,  gravity,  etc. 

THETAF-  Computes  potential  temperature  tendency  from  (12.11)  or  (12.15). 

THFCST-  Computes  the  potential  temperature  at  the  next  time  step  given  the 
tendency  from  THETAF.  The  forecast  values  are  space  smoothed. 
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TABLE  16.1  Start  Options  for  the  7L  ANSPE  Model 


DOCUMENTATION  OF  SETC  HOED 


FUNCTION 

"BITS" 

SETC 

BIT 

NUN 

STATE 

BIT  (1-36) 

13 

(0-35) 

12 

1 

Use  47*51  Fields  to  Initialise 

13 

12 

0 

Use  53*57  Fields  to  Initialise 

13 

13 

1 

Contingency  72  to  96HR  Forecast 

14 

13 

0 

No  72  to  96HR  Contingency  Forecast 

15 

14 

1 

Diagnostic  Print 

15 

14 

0 

No  Diagnostic  Print 

16 

15 

1 

Extended  48-72HR  PE  Forecast 

16 

15 

0 

Non  Extended  48-72HR  Forecast 

17 

16 

1 

Off-tine  PE  (06  and  18Z) 

17 

16 

0 

On-tine  PE  (00  and  12Z) 

18 

17 

1 

Update  (00-06) 

18 

17 

0 

No  update 

19,20,21,18, 

19,20 

3-BITS 

These  Bits  Are  Used  to  Conpute  the  Forecast 

22 

21 

1 

Start  Hour 

No  Internal  Start  of  PE12,  PE24... 

22 

21 

0 

Do  Internal  Start  of  PE12,  PB24,... 

23,24 

22,23 

0 

MAXTIM  -  48 

23,24 

22,23 

1 

These  bits  deternine 

MAXTIM  -  12 

23,24 

22,23 

2 

the  forecast  finish 

MAXTIM  ■  24 

23.24 

22,23 

3 

hour 

MAXTIM  -  36 

THFIX-  Time  smooths  the  potential  tempecature  and  geopotential  height. 

Moves  the  potential  temperature  (THETA  array)  and  geopotential  height 
(PHI  array)  forecasts  into  the  current  (N)  time  step.  Moves 
potential  temperature  and  geopotential  height  from  the  current  time 
step  into  the  last  time  step.  This  routine  is  called  by  UVFCST  to 
allow  u  and  v  at  the  future  time  step  to  be  stored  in  the  arrays 
THETA  and  PHI  after  the  pressure-gradient  force  is  computed. 

TIMSMP-  Time  smooths  the  P-Sigma  array  valid  at  time  N  to  avoid  time 

splitting.  The  array  valid  at  time  H  is  then  used  as  the  starting 
point  for  the  next  time  step;  that  is,  the  N  array  becoaws  the  N-l 
array  for  the  next  time  step. 

TIMSMT-  Same  as  TIMSMP  except  time  smoothes  the  potential  temper ature. 

TIMSMU-  Same  as  TIMSMP  for  the  u  and  v  wind  components. 

TIMSTP—  Uses  a  centered  time  step  to  compute  the  N+l  time  step  of  p  given  the 
tendency  value  and  value  at  the  N-l  time  step.  a 

TIMSTT-  Same  as  TIMSTP  for  potential  tender ature. 

TIMSTU-  Same  as  TIMSTP  for  the  u  and  v  wind  components. 

TNDTRK-  Truncates  the  tendencies  of  u,  v  and  potential  temperature  below  20 
degrees  north. 

UVFCST-  Computes  the  tendencies  of  u  and  v  (12.9)  and  (12.10)  or  (12.13)  and 
(12.14).  (Section  12.2} 

We  present  the  description  of  the  major  arrays  used  throughout  the  7L  AWSPE  model 

in  Table  16.2. 


Table  16.2.  List  of  Major  Arrays  Used  by  the  7L  awspe  Model. 


n  Blocks  FIELD1/FIELD2 


The  arrays  in  common  blocks  FXELD1  and  FIEL02  are  mirror  images  of  one  another. 

Common  block  FIELDl  contains  all  arrays  that  are  mapped  into  D-BANK1  and  common  block 
FIELD2  contains  all  arrays  that  are  mapped  into  D-BANK2.  The  array  names  in  FIELD2 
are  treated  as  dummy  names  and  are  never  referenced  in  the  model  equations,  but  the 
storage  locations  are  used  by  changing  the  bank  descriptor  index  (BDI ) . 


ARRAY  MAME 


DESCRIPTION 


CD  (53,  1,  1,  31) 

F  (53,  1,  1,  31) 
IMAX  (31) 

IMAXW  (31) 

IMINW  (31) 

JMAX 

JMAXM1 

JMAXM2 

JMAXW 

JMIN 

JMIN (31) 

JMINW 

JMINPl 

JMINP2 

JMNWP1 

JMNWP2 

JMXWM1 

JMXWM2 

PHI  (53,  7,  3,  31) 


PI  (53,  1,  3,  31) 
PIVAL  (3000) 

PSIG  (53,  2,  3,  31) 
RLNDSE  (53,  1,  1,  31) 
RM  (53,  1,  1,  31) 

RMAP  (53,  1,  1,  31) 


CD  (I,  1,  1,  J)  surface  drag  coefficient. 

F  (I,  1,  1,  J)  Coriolis  parameter. 

Maximum  I  index  for  half  grid  points,  function  of  31 
J  Values. 

Maximum  I  index  of  whole  grid  points,  function  of  31 
J  values. 

Minimum  I  index  for  whole  grid  points,  function  of  31 
J  values. 

Maximum  J  index  for  half  grid  points. 

JMAX-1 ' 

JMAX-2 

Maximum  J  index  for  whole  grid  points. 

Minimum  J  index  for  half  grid  points. 

Minimum  I  index  for  half  grid  points,  function  of  31- 
J  values. 

Minimum  J  index  for  whole  grid  points. 

JMIM+1 
JMIN+2 
JMINW+1 
JMINW+2 
JMAXW- 1 
JMAXW- 2 

PHI  (I,  K,  T,  J) ,  hgt  field  For  7  vertical  levels  and 
3  time  levels.  Also  used  to  hold  N+l  time  levels  of 
u  for  u  forecast. 

PI  (I,  1,  T,  J) ,  Exner's  function  for  three  time 
levels  and  one  vertical  level. 

PXVAL  (2*PRESSURE  (MB))  -  Value  of  PI  every  0.5  mb 
from  5  to  ISOOmb 

PSIG  (I,  KK,  T,  J) ,  p-sub-sigma  for  troposphere  and 
stratosphere  and  three  time  levels. 

RLNDSE  (I,  1,  1,  J) ,  sea  level  indcator  (0  for  sea, 
0.0001  for  land) 

RM  (I,  1,  1,  J) ,  map  factor  squared  with  xy  overbar 
applied  (half  grid  point  values) 

RMAP  (I,  1,  1,  J) ,  map  factor  squared  at  the  whole  . 
grid  points. 


SDOT 

(53,  6,  1, 

31) 

SHGT 

(53,  1,  1, 

31) 

THETA 

(53,  7,  3 

,  31) 

TSEA  (53, 

1. 

1,  31) 

U  (53,  7, 

2, 

31) 

P  (53,  1, 

1, 

31) 

V  (53,  7, 

2, 

31) 

VP  (53,  1 

,  1. 

31) 

ARRAY 

NAME 

WORK 

(53,  4,  1,  6 

WORK 2  (53,  2,  1,  61) 


ARRAY  NAME 


I HALF 


SOOT  (I,  K-l,  I,  J)  sigma  dot  for  levels  2-7, 
levels  1  and  8  are  zero. 

SHGT  (I,  1,  1,  J) ,  terrain  height  times  acceleration 
of  gravity. 

THETA  (I,  K,  T,  J) ,  potential  temperature  for  7 
vertical  levels  and  3  time  levels.  Also  used  to  hold 
N+l  time  level  of  v  for  in  UVFCST. 

TSEA  (1,  1,  1,  J) ,  sea  surface  temperature  (K) 
(climatological  values) . 

0  (I,  K,  T,  J) ,  u  time  levels  N-l  and  N. 

UP  (I,  1,  1,  J) ,  work  array  for  the  u  forecast. 

V  (1,  K,  T,  J) ,  v  at  time  levels  N-l. 

VP  (I,  1,  1,  J) ,  work  array  for  the  v  forecast. 

COMMON  BLOCK  FIELD 

DESCRIPTION 

WORK  (I,  IW,  1,  J+JJ) ,  where  IW  is  the  work  area, 

JJ* (IHALF-1) *30,  J«l,  31,  and  IHALF*2  for  dual  bank, 
version.  WORK  is  a  work  array  used  throughout  the 
model. 

WORK2  (I,  IW,  1,  J+JJ) ,  additional  work  area. 

COMMON  BLOCK  HALF 

DESCRIPTION 
DBANX  INDICATOR  (1  or  2) 


COMMON  BLOCK  PAR 


ARRAY  NAME 

DESCRIPTION 

BTHICK 

- 

Thickness  of  the  boundary  layer  (50  mb) . 

CP 

- 

Specific  heat  at  constant  pressure  (CP*1004.67) . 

CTDX 

- 

Delta  T/Delta  X  ■  2400  s/381000  m. 

DX 

- 

Delta  X-381000,M. 

G 

- 

Acceleration  of  gravity.  G*9.81  m/s. 

RKAPA 

- 

Constant  used  in  GEOPTL  (l./(l.+ROCP))*(1000)***ROCP 

R 

- 

Ideal  gas  constant  (R>287.04) 

RDLX 

- 

1/DX*1./381000 (1/ra) . 

RDLY 

- 

l/DY-l./381000(l/m) . 

RDT 

- 

1/DT-l. 2400 (1/8) . 

ROCP 

- 

R/CP. 

COMMON  BLOCK  IBARO 

ARRAY  NAME 

DESCRIPTION 

I  BAR 

- 

If  one,  run  barotropic  version. 

If  zero,  run  baroclinic  version. 
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COMMON  BLOCK  FIXER 


ARRAY  NAME 


DESCRIPTION 


KTMAX 

KTPMAX 

KPMAX 

IHMAX 


Maximum  number  of  vertical  layers. 

KTMAX+1,  maximum  number  of  vertical  levels. 

Maximum  number  of  P-sigma  layers. 

Maximum  number  of  data  banks  (the  maximum  allowable 
number  is  2) . 


COMMON  BLOCK  TONER 


ARRAY  NAME 


DESCRIPTION 


AFACL6 

APACG6 

RMUl 

TSMTH 


Pressure-gradient-averaging  factor  for  less  than  6  hours 
Pressure-gradient-averaging  factor  for  more  than  6  hours 
Smoothing  coefficient  for  latitudes  less  than  20N. 

Time  smoothing  coefficient. 


17.  SUMMARY 


In  Part  II,  we  have  discussed  the  numerical,  mathematical,  and  high 
level  program  design  o£  the  7L  AWSPE  model.  The  7L  AWSPE  is  a  dry  version  of 
the  NMC  7L  PE  described  by  Campana  (1978). 

He  developed  the  7L  AWSPE  for  two  reasons.  First,  operational 
requirements  dictated  moving  the  6L  AWSPE  from  System  V  to  System  I.  But,  the 
6L  AWSPE  could  not  meet  operational  timelines  when  running  on  System  I  because 
of  the  many  I/O  operations  required  by  the  model.  Second,  NMC  had  greatly 
improved  the  original  Shuman-Hovermale  model.  Campana's  work  made  it  possible 
for  AFGWC  to  run  a  fourth-order  version  of  the  7L  AWSPE  that  could  produce 
forecasts  almost  as  well  as  the  1/2-bedient  NMC  7L  PE  at  a  fraction  of  the 
cost. 


Like  Campana's  (1978)  version  of  the  NMC  7L  PE,  the  7L  AWSPE  uses  the 
hydrostatic,  sigma-coordinate,  energy-vorticity  form  of  the  primitive 
equations  on  a  polar  stereographic  projection  of  the  earth.  Two  finite- 
difference  formulations  of  the  model  exist;  a  second-order  and  a  fourth-order 
version.  The  model  domain  is  identical  to  the  6L  AWSPE  and  consists  of  a 
rectangular  area  centered  at  the  North  Pole.  At  the  outer  boundary,  a  no-slip 
boundary  condition  or  an  isothermal,  "nudge"  boundary  condition  can  be 
applied.  Terrain  is  included  in  the  model  through  the  sigma  coordinate 
system.  Time  and  space  smoothing  are  used  control  temporal  and  space  noise. 
The  only  thermodynamic  forcing  is  the  parameter isat ion  of  heat  flux  from  the 
ocean  surface  and  this  is  only  applied  at  the  lowest  layer.  Surface  drag  is 
parameterized  and  also  applied  to  the  lowest  layer  only.  Pressure-gradient 
averaging  is  used  to  allow  a  longer  time  step. 

The  7L  AWSPE  model  initialization  procedure  uses  the  nonlinear  balance 
equation  while  the  linear  balance  equation  is  used  by  the  6L  AWSPE.  The  7L 
AWSPE  model  forecasts  using  a  centered-time-dif ferencing  scheme.  The  output 
fields  are  identical  to  the  6L  AWSPE  model. 

The  7L  AWSPE  represents  a  transfer  of  technology  from  NMC  to  AFGWC.  When 
implemented,  the  fourth-order  version  of  the  7 L  AWSPE  model  will  produce 
significantly  improved  forecasts  of  height,  temperature,  and  wind  (compared  to 
the  second-order  6L  AWSPE) .  The  second-order  version  of  the  7L  AWSPE  should 
perform  only  slightly  better  than  the  6L  AWSPE.  Finally,  we  note  that  the  7L 
AWSPE  is  a  top-down  structured,  well  documented,  easily  maintainable  program, 
whereas  the  6L  AWSPE  is  none  of  these. 
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